 #
setClassUnion("numericOrNULL",members=c("numeric", "NULL"))
#class
setClass("FCMParamClusterAlg",contains="VParamClusterAlg",
		representation(
				.method="character",
				.dnoise="numericOrNULL",
				.alpha="numeric",
				.iter.max ="numeric",
				.nstart="numeric",
				.maxminJ="numeric",
				.m="numeric"
		),
		prototype=prototype(
				.description="FCM parametric Clustering algorithm class",
				.method="FCM",
				.dnoise=NULL,
				.alpha=0.001,
				.iter.max=100,
				.nstart=10,
				.maxminJ=5,
				.m=2
		
		),
		validity=function(.Object){
			if( !is.null(.Object@.dnoise) ){
				if(.Object@.dnoise <=0){
					stop("dnoise distance must be positive")
				}
			}
			if(.Object@.alpha >1){
				warning("alpha should be smaller than 1")
			}
			
			if(.Object@.iter.max <1){
				stop("iter.max must be greater than 1")
			}
			
			if(.Object@.nstart <1){
				stop("nstart must be greater than 1")
				
			}
			
			if(.Object@.maxminJ <1){
				stop("maxminJ must be greater than 1")
				
			}
			if(.Object@.m <1){
				stop("m must be greater than 1")
			}
		}
) 
#helper functions
#makes numeric cluster vector from cegclust output
#INPUT:	
#	clustObj - object generated by vegclust
makeClustNumVec <-function(clustObj){
	defuzzified <- vegclust::defuzzify(clustObj)
	return(as.numeric( as.factor( defuzzified$cluster ) ) )
}
clusterWithNoise <- function(fcmClustObject,inputSet,...){
	#TODO
	#take some parameters from fcmClustObject, estimate dnoise and call vclust function
	#at the end call clusterize function for new object
	
	#estimate initial noiseDist
	#TODO this should be square distances
	overallMeanVector <- colSums(inputSet$X)/dim(inputSet$X)[1]
	distances2MeanVector <- as.matrix( dist(rbind( overallMeanVector,inputSet$X) )^2 )[1,-1]
	noiseDist <- mean(distances2MeanVector)
	
	#estimate noise dist
	##initial values
	lambda <-1
	percentOfNoise <-0.0
	t <-1
	#TODO dnoise should be sqrt
	#print(paste("first RUN: ",fcmClustObject@.method))
	fcmClustsRes <- vegclust(x=inputSet$X, mobileCenters = fcmClustObject@.clusterNum 
			,method = fcmClustObject@.method, m = fcmClustObject@.m, dnoise = sqrt(lambda*noiseDist), 
			alpha = fcmClustObject@.alpha, iter.max = fcmClustObject@.iter.max,
			nstart=fcmClustObject@.nstart, maxminJ=fcmClustObject@.maxminJ)
	
	#min  exponent of double number
	VecsLen <- abs(.Machine$double.min.exp) +7
	noisePercentVec <- rep(0,VecsLen)
	noisePercentVec[1] <- calcNoisePerc(fcmClustsRes)
	lambdaVec <-rep(0,VecsLen)
	lambdaVec[1] <- lambda
	
	minDistBPoints <-min(dist(inputSet$X))
	
	while(percentOfNoise<=0.5){
		t <- t+1
		lambda <- lambda/2
		#print("Lambda");print(lambda)
		#extend lambda vec
		#lambdaVec <- c(lambdaVec,lambda)
		lambdaVec[t] <-lambda
		#calculate clusterization using centers frm previous clusterization
		#stopifnot(is.data.frame(fcmClustsRes$mobileCenters))
		#stopifnot( !is.na(sqrt(lambda*noiseDist)) )
	if(minDistBPoints > sqrt(lambda*noiseDist)){print("SMALL Delta!!!!!!!!!")}
		fcmClustsRes <- vegclust(x=inputSet$X,
				mobileCenters = fcmClustsRes$mobileCenters 
				,method = fcmClustObject@.method, m = fcmClustObject@.m, dnoise = sqrt(lambda*noiseDist), 
				alpha = fcmClustObject@.alpha, iter.max = fcmClustObject@.iter.max,
				nstart=fcmClustObject@.nstart, maxminJ=fcmClustObject@.maxminJ)
		
		#extend noisePercentVec
		percentOfNoise <- calcNoisePerc(fcmClustsRes)
		#noisePercentVec <- c(noisePercentVec, percentOfNoise)
		noisePercentVec[t] <- percentOfNoise
		
		#print("Percent of Noise");print(percentOfNoise)
		
	}
	
	#determine lambdaOpt using Pareto curve p ~ q*lambda^(-s)
	
	#find elements of noisePrecVec equal 0
	zeroElementsIndices <- which(noisePercentVec==0)
	if(!identical(zeroElementsIndices,integer(0))){
		#print("zero el")
		#print(zeroElementsIndices)
		noisePercentVec <- noisePercentVec[(-zeroElementsIndices)]
		lambdaVec <- lambdaVec[(-zeroElementsIndices)]
	}
	
	
	#after taking log we can use linear models fitting
	LOGnoisePercentVec	<- log(noisePercentVec)
	LOGlambdaVec		<- log(lambdaVec)
	
	
	#fit linear model 
	lambdaOpt<-1
	if(length(LOGlambdaVec)>1){
		lModel <- lm(LOGnoisePercentVec ~ LOGlambdaVec)
		s <- -lModel$coefficients[2]
		q <- exp(lModel$coefficients[1])
		lambdaOpt <- {q*s}^{1/{s+1}}
	}else{
		lambdaOpt <-1
	}
	#TODO S becomes negative
	if(is.na(lambdaOpt)){
# 		print(paste("Lambda:",lambdaVec,sep=" "))
# 		print(paste("Noise:",noisePercentVec,sep=" "))
# 		print(paste("S: ",s,sep=" "))
# 		print(paste("Q: ",q,sep=" "))
# 		print(paste("Method: ",fcmClustObject@.method,sep=" "))
		lambdaOpt <-1
	}
	#stopifnot(!is.na(lambdaOpt))
	fcmClustsResBest <- vegclust(x=inputSet$X,
			mobileCenters = fcmClustObject@.clusterNum 
			,method = fcmClustObject@.method, m = fcmClustObject@.m, dnoise = sqrt(lambdaOpt*noiseDist), 
			alpha = fcmClustObject@.alpha, iter.max = fcmClustObject@.iter.max,
			nstart=fcmClustObject@.nstart, maxminJ=fcmClustObject@.maxminJ)
	
	clusters <- makeClustNumVec(fcmClustsResBest)
	
	return(clusters)
	
	
}

#methods
#initialization method 
setMethod("initialize",
		signature="FCMParamClusterAlg",
		
		function(.Object,clusterNum,method="FCM",dnoise=NULL,alpha=0.001,iter.max=50,nstart=10,
				maxminJ=5,m=2,...){
			.Object <-callNextMethod(.Object,clusterNum,...)
			.Object@.method		<-method
			.Object@.dnoise		<-dnoise
			.Object@.alpha		<-alpha
			.Object@.iter.max	<-iter.max
			.Object@.nstart		<-nstart
			.Object@.maxminJ	<-maxminJ
			.Object@.m			<-m
			
			
			validObject(.Object)
			return(.Object)
		})


setMethod("clusterize",
		signature="FCMParamClusterAlg",
		definition=function(.Object,inputSet,...){
			#check input set
			callNextMethod(.Object,inputSet,...)
			print("Param FCM clust")
			clusters <-NULL
			if(.Object@.method == "NC" || .Object@.method=="GKN"){
				clusters <- clusterWithNoise(.Object, inputSet)
				
			}else{
				#print("vegClust start")
				#print(paste("Cluster NUM: \n", .Object@.clusterNum,sep=" "))
				clusterObj <- vegclust(inputSet$X,method=.Object@.method,mobileCenters=.Object@.clusterNum,
						m=.Object@.m,dnoise=.Object@.dnoise,alpha=.Object@.alpha,iter.max=.Object@.iter.max,
						nstart = .Object@.nstart,maxminJ = .Object@.maxminJ)
				clusters <- makeClustNumVec(clusterObj)
			}
			
			
			return(clusters)
		}
)
setMethod("get_description",
		signature="FCMParamClusterAlg",
		definition=function(.Object){
			paste(get_description(.Object),":",.Object@.method,sep=" ")
		})

